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A SIMPLE  MEANS  OF  UPDATING  THE  SRIF  FILTER  WHEN 
THE  STATE  EQUATIONS  ARE  IN  TRIANGULAR  FORM 

INTRODUCTION 

Estimating  the  state  of  a system  from  a set  of  uncertain  measurements  has  been  a prob- 
lem for  a long  time.  Kalman  in  the  early  sixties  provided  a simple  recursive  estimation  pro- 
cedure by  introducing  the  concept  of  state  and  state  transition.  This  procedure  in  some 
instances  provided  simpler  implementation  than  batching  techniques.  Since  Kalman’s  work  a 
number  of  numerical  procedures  have  been  developed.  An  excellent  account  of  these  pro- 
cedures as  well  as  historical  notes  can  be  found  in  Bierman’s  book  [1] . The  square-root 
information  filter  (SRIF  filter)  is  the  numerical  method  of  solving  the  Kalman-filter  equa- 
tions, v^ich  is  of  interest  in  this  report. 

There  are  a number  of  problems  which  involve  a state  transition  matrix  which  is  in 
upper  triangular  form.  Prominent  examples  of  problems  involving  the  condition  are  most 
tracking  problems.  This  report  describes  a simple  means  of  updating  the  prediction  process 
of  the  filter  under  this  condition.  A secondary  but  important  result  is  that  the  SRIF  filter 
lends  itself  to  parallel  hardware  implementation. 

REVIEW  OF  THE  SRIF  FILTER 

The  SRIF  filter  is  a numerical  method  of  implementing  the  Kalman  filter  [1] . The 
Kalman  fUter  is  obtained  from  modeling  the  process  as  state  equations,  defining  a measure- 
ment procedure,  and  best  estimating  the  states  of  the  systems.  The  state  equation  and  meas- 
urement process  are  defined  as 

X(k)  » ^ik)X(k  - 1)  + r{k)W(k) 

and 

XM(k)^H{k)X(k)+V(k), 

where  it  is  desired  to  best  estimate  Uie  n-by-1  state  vector  X{k).  The  remaining  quantities 
are  an  n-by-n  state  transition  matrix  *I>(fe),  an  n-by-p  matrix  r(fe),  an  m-by-n  measurement 
matrix  H{k),  and  an  m -by-1  measurement  Jtor  X^jik).  IF(fe)  and  V{k)  are  independent 
Gaussian  noises  with  the  properties 

E[wm-o, 

ElW(k)WU)]  =S(k)5j^, 

E[V(k)]‘0, 
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ElVik)V'U)]  =Q(k)8j^, 
and 

E[W{k)V'{j)]  = 0, 

where  6ji,  is  1 when;  ~ k and  is  0 otherwise.  The  covariance  matrices  S(k)  and  Q(k)  are  of 
dimension  p by  p and  m by  m respectively. 

The  best  estimate  of  X(k),  denoted  by  X(k)  in  the  standard  Kalman-filter  format,  is 


X(fe)  = X(k)  + K{knX„  (k)  -H{k)Xik)] , (1) 

where  K(k)  is  the  filter  gain,  given  by 

K{k)  = Pik)H'{k)Q-^{k),  (2) 

in  which  ^k)  is  the  smoothed  covariance  matrix,  with 

^^(k)  = p-^(k)+H’(k)Q-^{k)H(k).  (3) 

P(k)  is  the  predicted  covariance  matrix,  with 

P(k  + 1)  = 4>(fc  + l)P(fe)4)'(fe  + 1)  + r(Jfe  + l)S(ft  + l)r'(k  + 1),  (4) 

A 

and  X(k  1 ) is  the  prediction : 

A:(fe  + l)  = ‘I>(fe  + l)X(fe),  (6) 


The  filter  operates  in  a predict-and-correct  fashion.  This  suggests  a simple  derivation,  out- 
lined bdow. 

Equation  (1)  is  the  least-square  estimate  between  the  prediction  and  the  measurement 
at  the  kth  sample  which  is  obtained  by  minimizing  the  cost  function 

J(k)  - [X(fc)  - X(fc)]  'P-I(fc)tx(fc)  - X(k)]  + [Xj^ik) -HX(k)]  •Q-^{k)[XM{k)-HX{k)] 

(6) 

with  respect  to  X(fc).  The  value  of  X(k)  which  minimizes  </(fe)  is  denoted  by  X(k),  is  the  best 
estimate  of  X(k),  and  is  given  in  equations  (1)  throu^  (3).  Given  the  best  estimate  of  X(k), 
the  best  prediction  is  simply  equation  (5)  wifii  the  covariance  of  (4).  The  process  is  then 
simply  repeated  recursively,  with  equations  (4)  and  (5)  being  the  prediction  and  equations 
(1)  through  (3)  being  the  correction. 

Ihe  SRIF  filter  is  a means  of  implementing  the  Kalman  filter  which  depends  heavily  on 
Cholesky  decomposition  and  the  Householder  matrix  triangulation  algorithm  [1] . The 
Cholesky  decomposition  is  performed  on  a symmetric  positive-definite  matrix  by  factoring 
it  into  the  product  of  a lower  triangular  matrix  L and  its  transpose: 

Q-^LL' 

and 
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The  algorithm  for  obtaining  L,  found  in  reference  1,  is 


for;  = l, 

for  fe  =;■  + 1, n, 

and 


9ife  = ^ik  ~ tot  k=J  + l, ....  n and  I = k, ....  n. 

The  cost  function  in  equation  (6)  can  be  written  as 

J = (X-  X)'RR  '(X  -X)  + {Xm  - HXy(L')-^L-^  (X^  - HX),  (7) 

where  the  parenthetical  k has  been  dropped  for  notational  covenience,  P~^  is  factored  into 
RR',  and  Q(fe)  is  factored  into  L L'  (note  that  Q~^(k)  = (L')"^L"^ ).  Equation  (7)  can  be  re- 
written as 

J=(Z-R  'xy{Z  -R'X)*  {Zm  - H^XyiZ^  -H^X), 

where 


Z = R'X, 

Hyf  = L-^H. 


Equation  (7)  can  then  be  rewritten  more  compactly  as 


J = 


(8) 


The  cost  J is  unaltered  if  an  orthogonal  transform  T,  where  T'T  = /,  is  multiplied  by  the  new 
resulting  vector  in  (8).  Consequently  using 


= C 


1 

I 

1 


i 

1 

1 

1 

j 


I 

i 


inJ^&C  yields  the  same  cost  J as 


J=C'T'T  C. 
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In  addition,  if  T,  which  is  a n m square  matrix,  is  chosen  such  that 


and 


la  ■ ['.1 


then  the  cost  J becomes 


(9) 


(10) 


J = (R'X-Z)\R'X-Z)  + e'e. 

By  inspection  the  least-square  estimate  of  X is 

R'X^Z  or X = R-^Z, 

e'e  is  the  minimum  value  of  the  cost  J,  and  the  smoothed  covariance  is  P(k)  = (R')~^R~^ . 
For  simplification  (10)  is  argumented  to  (9),  yielding 


~R' 

1 

R'  Z 

Hw 

1 

0 e 

(11) 


The  transform  T triangulizes  the  matrix. 

To  show  (11)  is  equivalent  to  the  smoothing  portion  of  the  Kalman  filter,  (11)  can  be 
written  as 


1 

II 

'.R'  Z\ 

1 

Cl 

1 

1 

o 

1 

Equating  terms,  one  has 
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equation  (12)  becomes 


RR’  + HyfrH„=RR'. 

Using  previous  definitions  for  and  using P~^  =RR'  and  P~^  = R R',  one  obtains  equa- 
tion (3)  of  the  Kalman  filter.  Similarly,  substituting  (14)  and  (15)  into  (13),  one  obtains 
equation  (1)  of  the  Kalman  filter. 

The  Householder  algorithm  can  be  used  to  triangulize  the  matrix  represented  in  (11) 
without  ever  computing  the  transform  T directly.  Only  the  basic  results  are  sketched,  and 
an  example  is  given.  Detailed  information  may  be  found  in  reference  1.  The  algorithm  is 
based  on  reflection.  Let  the  vector  U be  normal  to  the  plane  U^.  An  arbitrary  vector  Y can 
be  represented  by 


Y=(Y'U)U*v,  (16) 

where  U=UI{U'  and  is  that  part  of  Y that  is  orthogonal  to  £7.  The  reflection  of  Y 
denoted  by  Y^  in  the  plane  is 

Y,  = -(Y'U)U  + u,  (17) 

and  the  results  are  represented  in  Fig.  1. 

Eliminating  v from  (16)  and  (17)  yields 

Yr=Y-2  U=(I-0UU')Y  = TY,  (18) 


where 


The  matrix  T is  an  elementary  Householder  transform  with  properties  T'  = T and  TT'  = I. 
Equation  (18)  can  be  shown  to  triangulize  a matrix  by  first  setting  the  elements  of  the  vec- 
tor Uby 
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u(l)  = y(l)  + a, 

«(2)  = y(2), 


and 


u(;)  = y(j), 


where  a = sgn  y(l)-^Y'Y.  The  transform  TY yields  y;.(l)  = a and  y^d)  = 0 for;  = 1,  2 

The  first  column  of  the  matrix  is  chosen  as  y(j)  in  order  to  set  u{j).  Equation  (18)  is  then 
applied  successively.  The  sign  on  (18)  is  changed  to  yield  positive  diagonal  elements  of  a, 
and  the  notation  P = 2/U'U  is  introduced.  The  algorithm  operating  on  successive  columns  of 
the  matrix  is 


Yr  = -Y  + p(Y'U)U. 


For  example 


«11 

“12 

' CO 

“1 

^2 

1 

eo 

«21 

“22 

“23 

0 

^22 

^23 

“31 

“32 

“33 

0 

^32 

^33 

“41 

“42 

“43 

0 

*>42 

*'43 

where 


U = 


Oji  +ai 


*21 


*31 


*41 


oi  =8gn(an)\/af^  +o|j  +a|i  +a^^, 

P = 2/U'U, 

Jlj  = fliy«(l)  + a2ju(2)  + 03yu(3)  + a^ju(4). 


and 


bjj  * - Oij  + PyijU(i)  for;  = 1,2,  and  3 and  i = 1, ...,  4. 


J 
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The  process  is  repeated  for  each  successive  submatrix.  The  next  step  is 


b\2 

*13 

“1 

*12 

*13 

1 

o' 

0 

*>22 

*23 

0 

“2 

“23 

0 

T2 

0 

*32 

*33 

0 

0 

“33 

• 

0 

*42 

*43 

0 

0 

“43 

where 

\*>22*^2 


'32 


02  = 8ign(b22)\A22  ■^*’32  ■^*’42* 
P = 2IU'U, 


72y  = b2jU(l)  + b2jU(2)  + b4yU(3), 


and 


Cjj  = - 6y  + 0rY2jU(j  - 1)  for/  = 2 and  3 and  i * 2,  3,  and  4. 


(19) 


Equation  (19)  is  the  desired  triangular  form  required  of  equation  (11)  for  the  example. 
The  correspondence  is 


* , 

“11  “12 

“1 

*12 

. R'^ 

“21  “22 

0 

“2 

“31  “32 

“41  “42 


r 

CO 

1 

Ci9 

“33 

, Z- 

, and  e = 

“23 

“23 

“43 

The  Householder  algorithm  just  described  can  be  compactly  encoded  in  Fortran  for 
general  computer  operation.  In  some  cases  a hardware  implementation  is  desirable  and  is 
shown  schematically  in  Fig.  2. 
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Fig.  2 — Schematic  of  operations  performed 
with  the  Householder  algorithm 


The  Kalman  and  SRIF  filters  were  briefly  reviewed  to  set  the  notation  and  acquaint 
those  readers  not  familiar  with  algorithms  with  the  salient  features.  A simple  means  of 
obtaining  the  prediction  portion  of  the  SRIF  Alter  under  an  important  special  case  is  next 
considered. 


PREDICTION  PROCESS 

The  smoothing  portion  of  the  Kalman  filter  using  SRIF  implementation  updates  the 
factorization  of  the  smoothed  covariance  and  the  transformed  best  estimate.  It  is  desirable 
to  update  the  prediction  process  in  a commensurable  form.  Only  an  important  special  case 
is  considered. 

The  process  noise  W{k)  is  assumed  to  be  zero,  and  the  state  transition  matrix  is 
assumed  to  be  in  upper  triangular  form.  Equation  (4)  updating  the  prediction  covariance 
then  becomes 
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where  the  noise  is  removed  and  the  sample  k has  been  dropped  for  notational  conveni- 
ence. The  inverse  of  (20)  is  taken,  yielding 

p-l  = (4)')-lp-l$-l  (21) 

The  covariances  are  replaced  with  their  factorization 

vdiich  can  be  rewritten  as 

Note  that  (^')~^R  is  in  lower  triangular  form,  which  means  that 

R » {^’)-^R  . (22) 

Equation  (22)  shows  the  simple  form  of  updating  the  factor  of  the  prediction  covariance. 
The  predicted  state  given  by 

X = ^X 

from  equation  (5)  is  transformed  by 

(R'y^Z  = ^(R'y^Z , 
where  X = (R'y^Z  and  X = (R'y^Z.  Solving  for  Z yidds 

Z’-(R')4>iR'y^Z.  (23) 

Substituting  R from  (22)  into  (23)  yidds 
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and  equation  (24)  remains  the  same  under  the  fading-memory  condition. 


IMPLEMENTATION 

As  an  example  a tracking  problem  is  taken  into  consideration.  The  state  transition 
matrix 


1 t 0 0 

1 -t  0 O' 

0 10  0 
0 0 1 t 

, where  = 

0 10  0 

0 0 1 -f 

0 0 0 1 

0 0 0 1 

represents  a target  moving  in  a straight  line  in  a two-dimensional  Cartesian  coordinate  sys- 
tem. The  components  of  the  state  vector  X(fc)  to  be  estimated  are  Xi  (k),  the  position  in  the 
ith  direction;  X2(I:)>  the  velocity  in  the  ith  direction;  X3(fc),  the  position  in  theyth  direc- 
tion; and  X^(k),  the  velocity  in  the ;th  direction.  Only  ^e  positions  are  measured;  conse- 
quently the  measurement  matrix  H is 


H = 


1 

0 


0 0 
0 1 


0 

0 


The  functional  flow  of  the  filter  is  shown  in  Fig.  3.  The  measurement  is  prewhitened 
using  the  Cholesky  factorization.  In  most  tracking  problems  the  inverse  required  can  simply 
be  written  in  closed  form  using  the  Cramer  rule.  The  prediction  variables  are  updated  wi^ 
no  more  than  a matrix  multiplication.  These  steps  can  be  mechanized  with  several  degrees  of 
parallelism  in  hardware.  Finally  the  smoothing  is  performed  using  the  Householder  algo- 
rithm shown  schematically  in  Fig.  2.  The  output  of  the  filter  in  normal  tracking  is  the  sta- 
tistical distance  [2,  3]  J = e'e  which  is  required  for  correlation  (a  direct  consequence  of  the 
filter)  and  the  predicted  position  X used  in  correlation  and  for  display.  The  outputs  are 
easUy  obtained,  induing  X,  because  need  not  be  found;^The  best  estimate  X can  be 
obtained  from  2 and  R'  directly  by  back  substitution,  since  iJ'  is  in  triangular  form.  All  the 
operations  described  including  the  Householder  algorithm  are  simple  operations  easUy 
mechanized  with  parallelism  in  the  hardware. 


SUMMARY 

The  SRIF  filter  was  briefly  reviewed,  including  the  Cbolesky  factorization  and  House- 
holder algorithm.  The  smoothing  portion  of  the  SRIF  filter  is  claimed  to  have  good  numeri- 
cal characteristics  and  lends  itself  to  parallel  hardware  operation.  The  prediction  process 
under  an  important  simple  case  was  examined.  The  state  transition  matrix  was  assumed  to 
be  in  upper  triangular  form,  and  the  process  noise  was  assumed  to  be  zero.  Most  tracking 
problems  can  be  formulated  in  this  form.  Under  this  special  case  it  was  shown  that  the  trans- 
formed smoothed  and  predicted  states  were  identical  and  that  the  smoothed  and  predicted 
covariance  factors  were  related  by  a simple  matrix  transform.  Consequently  the  entire  SRIF 
filter  including  both  the  smoothing  and  prediction  lends  itself  to  hardware  implementation. 
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Fig.  3 — Functional  flow  of  the  SRIF  filter 
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